Europhysics Letters 



PREPRINT 



Structure Below the Growing Surface 



E. Katzav 1 , S. F. Edwards 2 and M. Schwartz 3 

1 Department of Chemical Physics, The Weizmann Institute of Science Rehovot 76100, 
Israel 

2 Polymers and Colloids Group, Cavendish Laboratory, Cambridge University Mading- 
ley Road, CB30HE Cambridge, United Kingdom 

3 School of Physics and Astronomy, Raymond and Beverley Faculty of Exact Sciences, 
Tel Aviv University, Tel Aviv 69978, Israel 



PACS. 05.70.Ln - Nonequilibrium and irreversible thermodynamics. 
PACS. 02.50.-r - Probability theory, stochastic processes, and statistics. 
PACS. 81.10.Aj - Theory and models of crystal growth; physics of crystal growth, crystal 
morphology, and orientation. 



Abstract. - In recent years there has been a growing interest in the statistical properties 
of surfaces growing under deposition of material. Yet it is clear that a theory describing the 
evolution of a surface should at the same time describe the properties of the bulk buried 
underneath. Clearly, the structure of the bulk is relevant for many practical purposes, such 
as the transport of electric current in devices, transport of fluids in geological formations and 
stress transmission in granular systems. The present paper demonstrates explicitly how models 
describing deposition can provide us with information on the structure of the bulk. Comparison 
of an analytic model with a simulation of a discrete growth model reveals an interesting long 
range tail in the density-density correlation in the direction of deposition. 



The study of the statistical properties of surfaces growing under the deposition of material 
has attracted many researchers over the last two decades. The systems under consideration 
vary from heaps of sand or other assemblies of granular matter, to devices manufactured by 
the bombardment of atoms on a growing target. The theoretical description of such systems is 
given by a number of discrete and continuous models that belong mainly to three categories. 
The first is the Edwards Wilkinson category [1] which was constructed to describe a situation 
of slow deposition under gravity where each deposited particle has the time to find its lowest 
possible gravitational potential energy in the presence of the existing surface. The second 
category is that of KPZ [2] in which lateral growth is important. This can be a result of stick- 
ing or just the geometry of growth perpendicular to the surface as explained in ref. [2]. The 
third category is the MBE (Molecular Beam Epitaxy) [3] which was constructed to describe 
processes of device fabrication in which the physics should produce under a wide range of 
parameters flat surfaces. The focus of interest in those studies was the statistical character- 
ization of the growing surface. This is achieved by obtaining the roughness exponent of the 
steady state surface, the growth exponent [4-20] and the scaling functions associated with the 
steady state evolution of the surface [21-24]. It is easy to envisage many practical applications 
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for which the fluctuations in the steady state surface are relevant and this was always one of 
the motivations for the intensive study of surface growth. Yet it is obvious that the internal 
structure of the material below the surface may be of much more practical importance as this 
determines the mechanical properties of the system, generated by deposition, such as a heap 
of granular matter in which it will affect stress transmission. In electrical devices obtained by 
deposition, the structure of the bulk will obviously affect the transmission of electric current 
which will determine the functionality of the device. Certain geological formations relevant 
to the oil industry are also generated by deposition. The structure of the bulk determines the 
important flow properties through the formation. 

This letter has two goals. The first and more important goal is to draw attention to the 
fact that a theory which describes the evolution of the upper surface of a system growing 
under the deposition of material, should simultaneously be able to predict the structure of 
the bulk below that surface. The physical reason for our statement is the fact that serious 
rearrangement does not usually take place in the deeper layers below the surface. 

The second goal is to explicitly demonstrate how the structure of bulk below the surface 
can be obtained for two growth models, the continuous analytical KPZ system and the discrete 
numerical ballistic deposition (BD) system. Note, that the main reason for treating the two 
models mentioned above, is that both possesses as will be shown in the following share the 
common realistic physical property that the growth process is accompanied by the embedding 
of voids below the evolving surface (this property is not shared by other very successful discrete 
numerical systems such as SOS [25] and RSOS [5] which do not allow for the generation of 
voids below the surface.) The second reason for treating these two systems is that they are, 
widely known and simple enough to demonstrate the connection between the system describing 
the evolution of the surface and the structure of the bulk below it. 

To clarify our ideas consider first the Edwards Wilkinson equation for the local height of 
the surface h(r, t), 

^ = vV 2 h + R(r,t), (1) 

where R (r, t) is the local rate of deposition of material given by R (r, t) = Rq + r\ (r, t) and 
•q (r, t) is a noise term that has zero average and its correlations are usually taken to be 

(ry (r, t) rj (r', f )> =DA(r- r') 5 (t - t') , (2) 

where A is a short range function. The constant position and time independent rate i?o is tra- 
ditionally deleted, because it has no effect on the shape of the surface. For our purpose though, 
we have to keep it. To be precise note that the deposition rate is related to the local number 
of particles landing per unit area per unit time, n(r,t), by the relation R(r,t) = Qn(r,t) 
where D, is the effective volume taken by each landing particle. Now particles land and rear- 
range on the surface to minimize their potential (gravitational) energy. This rearrangement 
is described by the first term on the right hand side of eq. QJ. Dividing h(r,t) by fi we 
obtain the total area density of particles at the point r at time t. It is clear that whatever the 
dynamical picture of surface evolution is, the total number of particles has to be conserved 
and indeed eq. Q trivially conserves the number of particles. Actually it does more than 
that. It conserves the volume occupied by the particles. The picture is thus the following. 
The landing particles form a compact structure and then rearrange on the surface preserving 
the compact structure beneath it. If the particles are cubes, as in some discrete models, there 
remain no voids among the particles. Other shapes of particles must result, of course, in voids 
among them but the structure is expected to be either ordered or random close packing. Our 
following discussion is not concerned with that compact structure. We focus rather on devia- 
tions from that structure involving larger voids in the system. The physical reasons for such 
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Fig. 1 - The addition of the shaded particle creates a void in (a) and (c), where it sticks to the column 
on it right or left. No void is created (b). 



voids include sticking and the geometry of growth normal to the surface, that are described 
within the KPZ category. This incorporates voids into the structure in a most natural way. 
To visualize it, consider the following version of ballistic deposition. In this model a particle 
falls vertically and sticks to the first site along its trajectory that has an occupied nearest 
neighbor. The particle is not allowed to stick to a diagonal neighbor. The last particle to 
be deposited is shaded in Fig. JIJ that demonstrates how voids are created. A void, once 
created, does not disappear. 

In the following we will discuss the problem of structure below the evolving surface from 
two different points of view. The first is to show that within the continuous KPZ description 
of the evolving surface, it is possible to obtain the density-density correlation function of the 
material below the surface in terms of height correlations related to the evolution of the surface. 
The second is to obtain the density-density correlations directly from a one dimensional BD 
simulation. 

Consider next the KPZ equation 

dh 

j£ = vV 2 h + \(Vh) 2 + R(r,t). (3) 

The origin of the non linear term in the above equation is the fact that growth is perpendicular 
to the surface. It is clear that if a particle lands on the surface an outcrop perpendicular to the 
surface is generated. The overhang screens the area below it and prevents more descending 
particles to fill the void generated by that particle just as in the discrete BD version depicted 
in Fig. Q). The KPZ equation must be thus viewed not as an equation for a single valued 
height function, because the height function is not really single valued, but rather an equation 
for the height envelope function below which there is an abundance of voids. Can we see the 
existence of those voids from the KPZ equation itself? From the relation R(r,t) = fln(r,t) 
it is clear that the total " deposited volume" is given by J drdtR (r, t) (assuming we start 
deposition on a flat surface at time t = Q) . From eq. © on the other hand it is clear that the 
total volume below the envelop surface is larger than that by A J drdt (V/i (r, t)) which must 
be the volume occupied by voids. The average of the integrand is known in the literature as 
the excess velocity and its meaning is that the average height growth faster than expected 
from the rate of descending material and this is just due to voids being incorporated into the 
structure [11]. Suppose we could identify that part of the local rate of increase in the height, 

dhi Q t ^ , which is a result of void creation, then the local density p(r, z) at a point, through 
which the surface passed at time t, is given by 

dh (r,t) 



p(r,z) = po 



(dh(r,t) 


'dh (r,i)~ 




{ dt 


dt 


} 



at • < 4 > 
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where po is the constant density that would have existed in the corresponding Edwards Wilkin- 
son system. 

The left hand side of the above gives p as a function of z but on the right hand side we have 
functions of t. The next step must thus be to connect between the perpendicular coordinate 
z and the passage time of the surface t. This relation is readily obtained by noting that the 
height of the surface at the point r and time t is given by 



h(r,t) = R t + Sh (r, t) 



(5) 



where Sh(r, t) is the usual variable describing the width of the surface. This quantity can 
attain values of the order of a positive power in the lateral size of the system. Nevertheless, 
if we wait for long enough times, that make Rot ,to be of the order of the lateral size, the 
second term on the right hand side of eq. (JHJ is negligible compared to the first. The relation 
can be thus iterated in the following manner 



t — z/Ro — Sh (r, z/Ro — Sh(r,z/Ro •■■))■ 
Within the leading approximation to eq. © eq. is also simplified and reads 

~dh(r,t)~ 



We identify next 




dt 



Ro 



v,t—z/Ro I 



A(V/!(r,i)) 2 , 



(6) 



(7) 



(8) 



since this non-linear term is the term responsible for the "excess velocity" and thus for the 
incorporation of voids. Eqs. (JTJ) and JSJ will be used now to obtain the average density p and 
the density-density correlations in the system 



p (l-A((V^) 2 )/i?o) 



(9) 



Now, assuming that the steady state structure factor is given by (hqh-q) — Aq~ r for q < qo 
(where qo is the high momentum cut-off, corresponding to the size of landing particles) and 
zero above it, the average density is given by 



P = Po 



1 - XAS a 



(d+i-r) 



(d+l-T)R 



(10) 



where d+ 1 is the dimension of space and Sd is the surface area of a unit sphere in d dimensions. 
The density-density correlations involve one non-trivial correlation 



[Vh(r,t)} 2 [Vh(0,0)} 2 



(why 



(ii) 



To evaluate it we calculate it to the lowest order in the frequency dependent structure factor 
^(q, uj) — {hq^h-q,-^). Standard manipulation leads to the expression 



* ( r > *) = -~^2d A2 I dldV t 1 ' I ') 2 l~ T l'~ r f M) / M exp [i (1 + 1') 
(27r) J 



(12) 



where oji = Bl z is the typical frequency related to the decay of a disturbance of wave vector 
1 and the scaling function f(u) is related to the time dependent structure factor $ (q, t) = 
(27r) -1 / dw$(q,w) e luJt by 

$(q,t) = Aq- r f(Lo q t). (13) 
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Clearly, the calculation of higher orders in the iteration may be rather tedious but they are 
straightforward as described in refs. [13, 14]. From our point of view, however, the important 
statement is that density-density correlations below the growing surface are related to various 
height correlations which are natural objects of study in the traditional research of surface 
growth. It has to be carried in mind, though, that to use the correlation "J (r, t) for our 
purpose by relating the time difference to the difference in height, we must have Rot ^> 
\6h (r,0)\, 1^(0, t)|. 

We will work out as an example the one dimensional case and compare it with the sim- 
ulation of the one dimensional ballistic deposition. The density-density correlation function 
g (x, y) = (p (x) p (y)) — p 2 is related to ^ in the one dimensional case by 



(14) 



Using the scaling form for the ID time dependent structure factor [11], $ (q, t) = ^q~ 2 f (Bq 3 / 2 t) 
(where D is the noise amplitude, v the diffusion coefficient and B a numerical constant which 
cannot be determined analytically) we obtain 



8D 



1 



§tt 2 v 2 B a ' 3 t 4 / 3 



1 r / T x 2/3 

^73 /W cob U-) X 



dr 



(15) 



This form is based on the scaling form and is therefore correct for large t, so that for large y 
we obtain a power law behavior for 



«<»■'" = -(£ 



p8 



8D 



1 



9^2^2^4/3 



(y/Ro 



, 4 /3 



1 



.1/3 



/ (r) dr 



-4/3 



(16) 



It is interesting now to discuss briefly the case A < 0. The parameters Rq and A corre- 
sponding to our physical picture are positive, of course. The situation with positive Rq and 
negative A, which can be handled formally by our equations, has no easy natural interpreta- 
tion. Certainly this is not related to the RSOS models, which are supposed to correspond to 
negative A but produce compact structures. 

Let us turn now to the density-density correlations obtained from a one dimensional model 
of ballistic deposition. The model used is the nearest-neighbor ballistic deposition (NNBD) 
model, on a lattice with L = 1024. At each time step a column i is picked at random. A 
particle (just one particle!) falls vertically and sticks to the first site along its trajectory that 
has an occupied nearest neighbor as shown in Fig. ^ We found it simpler not to use periodic 
boundary conditions. Fig. [21 below presents the structure obtained after deposition. The 
alternating shades of grey correspond to particles deposited within different time intervals. 
The voids are the white regions. This figure corresponds to a lateral size of 512 sites. 

The long-range correlations obtained in the analytical calculation above can be already 
observed qualitatively in the plot of the cluster generated by the simulation (Fig. [2] below), 
where some the voids tend to extend upward for a long distance. In order to quantify this, 
we obtain the correlation g(0,y) from the simulation. We define the correlation function 



g{x,v) 



l^(p(r)p(r')) 



P 



(17) 



where x is in the lateral direction, y is in the perpendicular (growth) direction, r = (X, Y), 
r' = [X 1 , Y') such that \X -X'\=x and \Y-Y'\=y. The number of such pairs (r, r') is iV, 
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Fig. 2 - A BD cluster obtained by depositing 100, 000 particles on a substrate of size L — 512. A 
time step is denned by a deposition of a single particle. The different shadings correspond to different 
time intervals each corresponding to the deposition of 10, 000 particles. 



the average is over realizations of the randomness and p is the average density. Practically we 
take the average to be the average over runs. In this work we focused on correlations in the 
y (growth) direction since they exhibit nontrivial features such as algebraic tails. We depict 
g(Q,y) in Fig. |3| The total number of pairs that goes into the evaluation of g(0,y), which is 
the product of the number of pairs in one run times the number of runs is 2.4 x 10 10 . 

We depict the log of the absolute value correlation function, because the correlation func- 
tion drops fast such that it is difficult to observe certain features. For example, it is difficult 
to see that the correlations become negative. When we depict the logarithm of the absolute 
value we see the point where the function changes sign as a kink. (Because the y-axis is a 
discrete coordinate the correlations do not become zero but go from positive values where 
they are decreasing to negative value where they are increasing). The correlation drops very 
fast and may seem to be of short range. The analytical prediction though suggests that we 
should examine more carefully the tail of <?(0, y). In Fig. E^b) we depict a log-log plot of this 
tail, which reveals a power-law behavior or the tail, namely g (0, y) oc — y — 1 - 3Q1 . We see good 
agreement the between the numerical result given above and the analytical one. Both yield 
a negative tail decreasing in size as a power law with exponents that are quite close (1.3 vs. 
1.333). 




Fig. 3 - (a) The log of the correlation function in the growth direction, (b) A log-log plot of the tail. 
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We have chosen here to concentrate on the density-density correlation function, because 
it is the most obvious and traditional characterization of the structure. It is clear that for 
specific applications like stress transmission or transmission of electrical currents we might 
be interested in other attributes of the structure. The main point we make in this letter 
is that the same analytical and numerical techniques used for the study of the growth of 
surfaces can be used in the study of material structure below the surface. Once this point is 
realized other properties of the structure can be readily be obtained. In addition, we have 
obtained correlations in the growth direction with an unexpected algebraic tail. If this is not 
just an artifact of the two models studied, it may have far reaching practical and theoretical 
consequences. As a possible usage we claim that the anisotropy between the growth direction 
and the lateral growth is a signature of the growth direction. Thus, having at hand a specimen 
of a sedimentary rock, for example, our result suggest that one can determine its original 
growth direction by determining the principal direction at which an algebraic tail in the 
structure factor appears. 

We hope that the present work will trigger more research in this direction, both experi- 
mentally, and by introduction of more realistic discrete models of deposition. 

* * * 
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